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1 3.  ABSTRACT  (Minimum  20 0  word t) 


In  this  report,  the  atmospheric  conditions  are  found  by  knowing  only  the  projectile's  flight  trajectory  and  its  flight 
coefficients  together  with  the  initial  atmospheric  conditions  on  the  ground.  The  test  trajectories  were  generated  as 
solutions  of  the  modified  point  mass  (MPM)  equations  of  motion.  The  correct  atmospheric  conditions  for  the  generated 
flight  trajectory  arc  obtained  from  data  collected  during  a  weather  balloon  flight.  A  nonlinear  least  squares  method  was 
then  used  to  lit  the  MPM  equations  to  the  lest  trajectory  by  varying  the  meteorological  parameters.  Density,  sound  speed, 
and  wind  profiles  agreed  well.  Further  tests  of  the  method  involved  the  flight  coefficients,  which  we  perturbed  to  observe 
the  corresponding  variation  for  the  fitted  values.  Also,  Gaussian  noise  was  introduced  onto  the  trajectory  values  to 
simulate  the  uncertainty  in  trajectory  measurements.  The  results  of  the  analysis  shows  that  the  whole  trajectory  should  be 
fitted,  and  not  small  trajectory  segments,  to  obtain  accurate  atmospheric  parameters  with  the  location  precision  provided 
by  current  radar  and  Global  Positioning  Satellite  System  (GPS)  techniques. 
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1.  INTRODUCTION 


A  better  knowledge  of  atmospheric  conditions  could  improve  the  accuracy  of  artillery. 
Currently,  weather  balloons  are  used  to  gather  the  temperature,  pressure,  and  wind  velocity 
as  a  function  of  height  above  the  earth’s  surfa-e.  As  the  weather  changes,  the  information 
becomes  stale  and  contributes  significantly  to  the  artillery  round’s  error  budget  (Matts  and 
Ellis  1990).  IdeaHy,  the  atmospheric  conditions  need  to  be  known  immediately  before  a 
round  is  fired. 

Since  the  trajectory  of  a  projectile  depends  upon  the  atmospheric  conditions,  among 
other  things,  the  problem  can  be  turned  around  to  determine  the  atmospheric  conditions 
by  knowing  the  trajectory  of  a  fired  projectile.  More  specifically,  a  nonlinear  least  squares 
method  can  be  used  with  the  trajectory  data,  if  the  equations  of  motion  are  assumed  with  a 
knowledge  of  the  projectile’s  flight  coefficients,  to  yield  a  best  fit  to  the  density,  temperature, 
and  the  wind  velocity  as  a  function  of  altitude.  Thus,  ♦he  atmospheric  conditions  would  be 
determined  soon  after  the  first  shot  was  fired,  the  elapsed  time  for  the  extraction  of  the 
meteorological  results  depending  on  the  computer’s  speed  and  the  efficiency  of  the  solution 
algorithm.  The  atmospheric  quantities  could  subsequently  be  used  with  the  artillery  piece 
swung  to  a  new  azimuth. 

As  a  first  step,  Cooper  and  Bradley  (1991)  posed  a  problem  that  did  not  use  atmospheric 
data  from  a  balloon  flight;  to  use  complicated  experimental  data  would  pose  a  harder  prob¬ 
lem.  Their  “observed”  trajectory  was  calculated  from  the  modified  point  mass  (MPM) 
equations  of  motion  (Lieske  and  Reiter  I960,  Bradley  1990)  combined  with  the  1966  U.  S. 
Standard  Atmosphere  (STAT)  and  linear  wind  profiles.  The  STAT  assumes  a  linear  temper¬ 
ature  with  altitude  and  a  simple  relationship  between  the  temperature  and  density.  To  facil¬ 
itate  fitting  the  manufactured  trajectory  data,  Bradley  (1990)  revised  the  MPM  equations 
so  that  current  fitting  techniques  could  be  adopted  to  iterate  the  air  temperature  gradient 
with  respect  to  height  and  the  wind  velocity  as  tittable  parameters.  FIN  LIE  (Bradley  1981), 
a  nonlinear  curve-fitting  method,  then  fitted  che  revised  MPM  equations  to  the  generated 
data  by  adjusting  these  parameter  values.  When  the  fitted  values  for  trajectory,  air  density, 
air  temperature,  and  wind  profiles  were  compared  with  the  generated  corresponding  values, 
the  agreement  was  found  to  be  excellent. 

The  program  to  implement  the  least-squares  fit  is  quite  cornplic;  ted  and  uses  a  large 
amount  of  time  on  a  computer’s  central  processing  unit  (CPU).  A  FORTRAN  subroutine  was 
constructed  for  F1NLIE  that  defines  the  original  set  of  equations  plus  an  auxiliary  set.  Tire 
auxiliary  equations  involve  partial  derivatives  of  the  parameters  to  be  fitted.  Some  of  these 
equations  are  exceedingly  complicated  and  the  chances  of  error  in  obtaining  these  equations 
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are  large  if  paper  and  pencil  were  to  be  used.  To  minimize  errors,  a  software  package  called 
MACSYMA  (MACSYMA  1983)  was  used.  MACSYMA  determined  the  partial  differential 
expres  ions  from  the  original  equations  and  translated  them  to  readily  usable  FORTRAN 
code. 

In  another  approach  that  has  not  yet  been  implemented,  the  temperature  and  density 
from  the  most  recent  atmospheric  data  are  assumed  correct  and  the  shell  trajectory  is  then 
used  to  calculate  the  wind  velocities  (Pedersen  1994).  Since  only  the  wind  velocity  needs 
to  be  determined,  the  calculations  should  proceed  rapidly.  The  possibility  exists  that  the 
temperature  and  density  could  have  been  measured  a  few  hours  ago,  thus  introducing  added 
uncertainties  with  time.  Possible  ways  of  modeling  the  density  and  temperature,  given 
changes  in  the  ground  meteorological  values,  are  being  explored. 

The  work  described  in  this  report  follow  on  and  expand  that  of  Cooper  and  Bradley 
(1991).  Instead  of  using  an  idealized  atmospheric  model  with  the  wind  velocity  being  linear 
with  height  above  the  ground,  the  atmospheric  data  from  a  weather  balloon  were  used  as 
input  to  generate  the  test  trajectory.  The  test  trajectory  was  then  fitted  over  segments 
to  obtain  approximations  to  the  wind  velocities,  density,  and  temperature.  Error  studies 
were  also  performed.  In  practice,  the  flight  coefTi  >ents  H'e  1  mown  only  approximately.  This 
uncertainty  was  studied  by  perturbing  the  accepted  flight  coefficient  values  to  see  how  the 
meteorological  values  varied.  Also,  white  noise  was  introduced  onto  the  trajectory  values  to 
simulate  the  uncertainty  in  trajectory  measured  by  radar  or  the  global  positioning  satellite 
system  (GPS)  techniques. 

2.  TRAJECTORY  EQUATIONS  AND  COORDINATE  SYSTEM 

The  construction  of  firing  tables  depends  significantly  upon  the  use  of  the  MPM  trajectory 
model.  More  complex  models  are  used  only  for  special  cases,  and  the  complete  data  ret  to 
construct  these  more  complex  models  is  not  commonly  available.  The  6-degree-of-freedoin 
equations  have  time  derivatives  that  do  not  appear  in  an  explicit  factored  form.  FIN  LIE, 
however,  uses  only  equations  with  derivatives  that  are  factored.  Thus,  the  MPM  equations 
will  be  used  here  in  a  different  but  almost  equivalent  form  (Bradley  1990)  than  is  commonly 
used  in  the  Firing  Tables  Branch  (Lieske  and  Reiter  1900). 

Bradley’s  formulation  for  the  MPM  model  can  be  re-written  in  the  form 

0  =  (D-  ,l)V  +  (7/t  (1) 
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and  in  which  U  —  projectile  velocity  with  respect  to  the  earth 
V  =  U  —  W  =  projectile  velocity  with  respect  to  the  air 
W  =  wind  velocity  with  respect  to  the  earth 
G  =  the  sum  of  the  gravity  and  Coriolis  accelerations 
<f>  —  axial  spin,  rad/s. 

(All  symbols  are  defined  in  the  List  of  Symbols.)  The  axial  spin  <f>  is  obtained  as  the  solution 
of  a  simplified  roll  equation: 

*  =  -BCt,  (£)  (2) 

in  which 


The  aerodynamic  coefficients  Cx0)  Cm0,  (wvpa,  and  Ct  in  Eqs.(l)  and  (2)  are  tabulated 
as  functions  of  Mach  number.  The  drag  coefficient  Co  depends  on  both  Mach  number  and 
the  yaw  of  repose,  a*e.  In  particular,  Co  is  assumed  to  have  the  form 


Co  =  Coo  +  Co2|oe|2 


(3) 


where  Coo  is  a  function  of  Mach  number  and  Cori  is  a  constant. 


The  yaw  of  repose  can  be  computed  from  «he  relation 
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A  convenient  coordinate  system  is  needed  for  describing  the  motion  of  a  projectile  along 
its  trajectory.  Following  Cooper  and  Bradley’s  (1991)  treatment,  we  assume  that  the  launch 
point  is  at  sea  level.  Then  we  set  our  origin  at  the  launch  point.  We  then  define  a  right- 
handed  Cartesian  system  as  follows:  the  1-  and  3-axes  form  a  plane  tangent  to  the  earth  at 
the  origin;  the  2-axis  is  perpendicular  to  this  plane,  positive  upward,  and  the  l-axis  is  chosen 
so  that  the  velocity  U  at  time  zero  is  in  the  1-2  plane.  Then,  the  projectile’s  position  vector 
X  with  respect  to  the  earth  can  be  written  in  component  form  as 


X  =  (XU  x2,  *3) 


in  which  X\  is  the  down-range  distance,  X2  is  the  height  above  the  1-3  plane,  X3  is  the  lateral 
distance,  positive  to  the  right  when  looking  down  range,  and  X  —  U.  Since  the  trajectory 
(in  most  cases)  lies  nearly  in  the  1-2  plane,  X3  is  usually  much  smaller  in  magnitude  than 
X\  and  X2  • 


Similarly,  we  can  write 

U  =  (Uu  t/2,  U3) 

V  =  (Vt,  V2,  V3) 

W  =  (Wu  W2,  W3) 

G  —  (G 1,  G2,  G3). 

in  which  U3  is  usually  much  smaller  in  magnitude  than  U\  and  U2.  The  initial  velocity  is 
given  by 

U0  =  1 1 ( cos  E ,  sin£,  0) 

in  which  E  is  the  gun  elevation.  The  details  of  the  launching  point  are  covered  in  the  report 
of  Cooper  and  Bradley  (1991). 


Eq.(l)  can  be  written  in  component  form  as 
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(5) 

(6) 
(7) 


The  value  of  G  is  the  sum  of  the  gravity  and  Coriolis  accelerations: 


G  =  g  +  C 


(8) 
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The  solution  scheme  used  for  this  study  has  one  important  modification.  Formerly,  tin- 
influence  equations  did  not  include  the  yaw  of  repose;  for  this  study,  the  influence  equations 
must  include  the  yaw  of  repose. 


3.  ATMOSPHERIC  DATA  AND  MODELING  OF  ATMOSPHERE  FOR 
FITTING  DATA 

In  this  study,  we  use  the  atmospheric  conditions  that  are  obtained  from  a  weather  balloon. 
We  hope  that  the  atmospheric  model,  which  is  assumed  to  fit  the  data,  is  close  enough  to 
reality  to  allow  good  approximations  of  the  measured  conditions.  More  complex  models  take 
more  computer  time  and  sometimes  do  not  yield  as  good  approximations  as  the  more  simple 
model.  The  latitude  L  for  the  simulated  test  was  39.15°  North,  the  azimuth  AZ  was  21.915°, 
and  g0  was  9.80665 m/s2. 

3.1  Atmosphere  U«ed  in  the  Generating  Equations  To  ge.ner  te  test  trr  Recto¬ 
ries,  we  used  the  balloon  data  taken  at  Aberdeen  Proving  Ground  (APG)  on  30  March  1992 
in  the  early  afternoon.  The  data  consisted  of  pressure,  temperature,  and  the  wind  velocity 
in  the  horizontal  plane  taken  at  very  small  intervals  of  altitude  initially  and  increasing  to 
intervals  of  400  meters  above  altitudes  of  3000  meters.  The  wind  velocity  data  are  shown  in 
a  later  section  when  comparisons  are  be  made  with  fitted  values. 

The  temperature  as  a  function  of  altitude  is  given  in  Figure  1. 


Figure  1.  Temperature  Data  Compared  to  Standard  Atmosphere  Curve  Temperature 


Shown  for  comparison  is  the  temperature  curve  obtained  from  the  U.S.  Standard  Atmosphere, 
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1966 ,  hereafter  railed  STAT  The  SI  AT  rurve  ran  he  nussub-red  as  an  vipproximat  n<ti  of 
the  average  temperature  in  temperate  inmates  in  N- »ri h  \meri<a  Although  the  recorded 
temperature  varies  no  more  than  a  few  dcgn*es  from  the  S  1  A  I  values.  the  temperat me  prottle 
can  vary  markedly  according  to  location  and  t  line  of  the  year.  I  he  winter  temperat  ure  profile 
in  Siberia  would  have  an  initial  positive  gradient  with  respect  to  alt  it  tide,  whereas  t  he  average 
temperature  profile  at  the  latitudes  in  North  America  would  have  a  temperature  gradient  of 
opposite  sign. 

The  density  data  arc  shown  in  Figure  2,  together  with  the  density  curve  obtained  from 
STAT  profile  for  locations  at  latitudes  in  North  America.  The  density  data,  although  not 
measured  directly,  are  calculated  from  the  pressure  and  temperature  data.  The  trial  input 
trajectory  was  obtained  by  numerically  integrating  the  MPM  equations  with  the  interpolated 
balloon  data. 


Figure  2.  Density  Data  Compared  to  Standard  Atmosphere  Curve  Density 


3.2  Atmosphere  Mode!  used  with  the  Fitting  Equations  Cooper  and  Bradley 
(1991)  fitted  the  horizontal  components  of  the  wind,  as  well  as  the  density  and  temperature, 
in  one  segment  with  a  linear  curve.  Not  surprisingly,  they  obtained  a  good  fit  with  the 
imposed  linear  wind  velocities.  For  the  present  work,  we  fitted  the  quantities  within  a  few 
time  segments.  Initially,  we  fitted  a  linear  relationship  for  the  wind  velocities  within  each 
segment,  but  we  found  that  parameters,  Cy,  of  constant  value  throughout  each  time  segment 
gave  a  sufficient  fit  with  the  advantage  of  using  fewer  fitting  parameters.  Essentially,  for 
each  segment  along  the  trajectory,  we  have  for  the  wind  velocity, 
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H . 

=  r, 

(») 

ir2 

-  0 

(10) 

H', 

=  r2 

(ii) 

The  rate  of  change  of  temperature  with  respect  to  altitude,  //,  is  the  third  parameter 
fitted  for  most  of  this  study,  C3  =  T' .  Here,  T'  is  the  lapse  rate  with  opposite  sign.  Cooper 
and  Bradley  (1991)  assumed  a  temperature  dependence  on  the  height  of  the  form 


T  =  Tb  +  T\H-Hb)  (deg  K)  (12) 

in  which  the  gradient  of  the  temperature,  T\  is  the  negative  of  the  lapse  rate.  The  value 
H  is  called  the  geopotential  altitude  and  takes  the  rotation  of  the  earth  into  consideration. 
The  geopotential  altitude  is  discussed  in  more  detail  by  Cooper  and  Bradley  (1991).  The 
value  Hb  is  the  value  of  //  at  the  bottom  of  a  geopotential  zone  for  the  STAT.  For  those 
zones  in  which  T'  is  not  zero,  the  STAT  air  density  is  given  by 


JL 

PD 
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T  yl-Q/T' 
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in  which  the  subscript  B  denotes  the  value  of  a  variable  at  H  =  Hb ,  and 

9o  M 


Q  = 


n 


(13) 


in  which  Tt  is  the  universal  gas  constant  and  M  is  the  molecular  weight.  The  value  of  Q  for 
dry  air  is 


Qi  =  34.16319474  (deg  K/km) 

For  moist  air,  the  value  of  Q  can  be  approximated  as 

Q  — - Qi - 

1+0.61/? 

The  value  /?  is  the  mass  fraction  of  water  vapor  in  the  air.  This  value  of  /?  can  be  as  high 
as  0.04.  For  this  study,  we  assumed  that  /?  =  0.  A  different  expression  is  used  in  which  the 
lapse  rate  is  zero,  but  for  this  study,  the  lapse  rate  is  never  zero. 


For  the  current  study,  the  lapse  rate  and  wind  velocities  were  fitted  on  altitude  zones 
assuming  that  the  density  varied  as  Equation  (13).  The  subscript  values  at  the  bottom  of 
the  corresponding  altitude  zones  are  denoted  as  A.  As  before,  the  first  time  interval  values 
of  T4,  H, 4,  and  pa  (obtained  by  measurements  at  the  launch  site)  are  required  inputs  to  the 
fitting  process;  thereafter,  the  closing  values  for  the  k-th  interval  become  the  starting  values 
for  the  (k+l)-th  interval. 
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4.  CHARACTERISTICS  OF  THE  TEST  PROJECTILE 


We  used  an  M107  projectile  to  generate  trajectories.  The  M107  projectile  has  been  well 
studied  and  all  the  coefficients  needed  for  use  with  the  MPM  equations  are  well  known 
(MacAllister  and  Krial  '975).  This  projectile  has  the  following  physical  properties: 

(  (diameter)  =  155  mm 

m  (mass)  =  43.09  kg 

I x  ~  0.1461  kg  —  m2 

Of  the  six  aerodynamic  coefficients  involved  in  our  equations  of  motion, 

Cdo,  C'o2,  Cia,  C\fa,  Ch’pa ,  and  C(p, 
four  are  functions  of  Mach  number,  as  indicated  in  Table  1. 


Table  1.  Mach  Dependency  of  the  Test  Case  Aerodynamic  Coefficients 


Mach 

Cdo 

Mach 

cLa 

Mach 

CjWo 

Mach 

Ct, 

0.750 

0.1  i50 

0.6 

1  650 

0.69 

3.30 

0.50 

-0.0145 

0.825 

0.1200 

1.9 

1.834 

0.74 

3.33 

1.00 

-0.0123 

0.875 

0.1310 

1.5 

2.294 

0.79 

3.40 

2.00 

-0.0096 

0.910 

0.1545 

3.0 

2.527 

0.85 

3.75 

3.00 

-0.0079 

1.025 

0.3800 

0.93 

4.67 

1.050 

0.4170 

0.97 

3.90 

1.100 

0.3900 

0.99 

3.70 

1.250 

0.3660 

1.04 

3.56 

1.600 

0.3250 

1.10 

3.51 

2.000 

0.2900 

3.00 

2.61 

3.000 

0.2200 

This  table  shows  pairs  of  values:  Mach  number  and  corresponding  coefficient  value.  Our 
code  performs  straight-line  interpolation  for  Mach  values  between  two  entries.  The  remaining 
two  aerodynamic  coefficients  are  assumed  to  be  constant: 

Cd2  =  4.C 
CNpa  =  -.75 
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The  muzzle  velocity  for  the  M107  projectile  was  397.4  m/s.  The  initial  spin  was  taken 
to  bo  H48  rad/s. 

Exactly  the  same  aerodynamic  behavior  was  assumed  in  the  fitting  equations  as  was  used 
in  the  equations  for  generating  the  trajectories.  For  this  caliber  shell,  it  is  generally  conceded 
that  the  drag  coefficient  is  known  to  within  a  half  percent.  The  significance  of  such  an  error 
will  be  explored  in  the  simulations. 


5.  PRELIMINARY  ANALYSIS 

The  !cr  t  squares  method  needs  a  unique  set  of  parameters  that  will  yield  a  global 
minimum.  A  first  consideration  is  whether  the  parameter  solution  of  the  flight  dynamic 
equations  are  unique  The  value  A  in  Equation  (1)  contains  Cp  and  the  parameter  p  as 
cofactors.  The  value  of  Cp  depends  on  the  temperature,  which  has  a  parameter  to  be 
determined.  Therefore,  the  value  A  would  seem  to  contain  two  factors  that  can  vary  and 
still  leave  the  value  of  A  unchanged.  However,  the  value  of  T  depends  upon  p  through  a 
constraining  equation,  and  uniqueness  should  be  restored.  Also,  the  parameter  Cp  is  weakly 
influenced  by  other  variables  since  Cp  depends  upon  the  yaw  of  repose. 

Another  consideration  is  the  sensitivity  of  the  parameter  values  to  errors  in  measurement. 
Consider  the  shell  at  apogee,  where  the  shell  velocity  is  aligned  approximately  along  the  X\ 
direction.  With  these  conditions  and  using  the  notation  of  Cooper  and  Bradley  (1991),  the 
MPM  model  yields  much  simplified  expressions  which  are,  in  component  form, 


Ux  = 

-AV\ 

(14) 

u2  = 

G2 

1  +  ha 

(15) 

u3  = 

HlGi 

(1  +  ^o)(l  —  h\f) 

(16) 

Assuming  flat  fire,  v’e  transform  from  the  time  to  the  X\  coordinate  and  focus  our 
attention  on  Equation  (14) 

U[  =  -a(l-W1/Ul)(Ul-Wl)  (17) 

in  which  a  —  A/(V/£).  This  equation  can  be  further  simplified  with  the  assumptions  of  no 
winds.  Transforming  to  nondimensionalized  variables,  x  =  (Xi  —  X\b)/£  and  i  =  U\b(t  —  tb)/£, 
and  integrating  we  obtain 

i=lnO±«0  (,8) 
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(Jib  is  the  value  of  the  first  component  of  the  velocity  at  the  beginning  of  the  trajectory 
segment  and  tb  is  the  value  of  t  at  the  beginning  of  the  segment.  For  a  small  enough  segment 
the  above  equation  can  be  approximated  by  expanding  in  series  to  two  terms, 

x  =  <(1  -  ai/2).  (19) 

The  question  of  the  accuracy  with  which  one  can  determine  p  can  be  examined  by  using 
Eq.  (19).  Many  of  the  measurement  techniques  smooth  the  output  positions  and  velocities 
so  that  one  is  not  aware  of  the  fundamental  error  spread  of  the  basic  data.  Given  that  this 
error  in  starting  and  ending  positions  is  Ax,  we  want  to  find  the  required  distance,  x,  of  the 
segment  to  make  the  relative  density  error  smaller  than  a  certain  value,  or  more  generally, 
the  relative  error  in  a.  We  differentiate  the  logarithm  of  Equation  (19)  to  obtain 

Ax  P  Aa 

~  =  _  2(t  —  aP/2) 

The  quantity  in  parantheses  in  the  denominator  on  the  right-hand  side  is  approximately 
equal  to  x  for  sufficiently  small  t.  Perform  the  substitution  and  rearrange  to  obtain 

x  =  \J\1.  Ax/  a  (Aa/a)|.  (21) 

Here,  (A a)/a  is  the  relative  uncertainty  in  a  that  can  be  tolerated,  and  Ax  is  the  uncertainty 
in  x  attributable  to  measuring  system  limitations. 

The  value  of  Cd0  is  generally  not  known  to  be  better  than  a  half  percent.  This  uncertainty 
is  partly  attributable  to  limitations  of  measuring  equipment  for  the  aerodynamics  range  and 
partly  to  round-to-round  differences.  Thus,  the  relative  error  of  a  cannot  be  less  than  a  half 
percent.  The  value  of  x  is  a  minimum  value  that  v/ill  give  A  a/ a  the  allowed  precision.  A 
representative  value  of  a  for  the  M107  projectile  at  ground  level  for  the  numerical  experiment 
is  a0  =  1.3  •  10"5.  The  drag  coefficient  decreases  as  the  shell  ascends  and  decelerates. 
Likewise,  the  density  decreases  with  altitude  and  the  value  of  a  decreases  accordingly.  The 
resulting  envelope  for  the  minimum  uncertainty  for  the  trajectory  length  in  which  a  relative 
uncertainty  of  a  is  <  1%  is  shown  in  Figure  3.  For  the  test  trajectory  at  or  near  apogee,  the 
value  of  a  is  only  about  a  fifth  of  its  value  on  the  ground.  In  practice,  the  relative  uncertainty 
in  a  will  be  on  the  order  of  the  relative  uncertainty  in  the  quantity  p ,  which  is  to  be  found. 

From  Figure  3,  it  can  be  seen  that  the  segment  lengths  for  the  M107  ( i  =  0.155 m)  must 
be  x  >  13,000  (13,000  x  0.155  m  (=»  2000  m)  if  the  uncertainty  in  the  segment  length  is  1 
meter  (A  x  =  1  m/0.155  m)  near  apogee  where  the  density  times  the  drag  coefficient  is  only  a 
fifth  of  what  it  is  at  launch  ( a/aQ  =  0.2).  The  length  of  such  a  segment  is  large,  and  the  shell 
may  change  altitude  more  than  is  desired.  To  achieve  this  much  precision  in  measuring  the 
position  may  be  a  daunting  task.  Equation  (21)  is,  of  course,  relevant  to  the  spacing  of  the 
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data  stations  in  a  firing  range.  For  a  7.62-mm  projectile  with  the  flight  characteristics  of  an 
M107  shell,  the  segment  length  would  need  to  be  approximately  1  meter  to  obtain  a  relative 
accuracy  of  0.005  in  the  drag  coefficient.  The  range  stations  in  the  ARL  Aerodynamic  Range 
are  placed  at  least  a  meter  apart.  Moreover,  th?re  are  many  stations  along  a  considerable 
length.  These  features  allow  the  Aerodynamics  Range  to  8chieve  the  desired  accuracy. 


(Aa)/a  =  0.01 


25000 

20000 

15000 

10000 

5000 

0.2 


a/aQo .  6 


Figure  3.  Required  Segment  length  (cal)  to  Achieve  Given  Accuracy  in  a  as  a  Function  of 
a  and  Uncertainty  in  Segment  Length 

The  relative  effects  of  the  density  versus  the  effects  of  the  change  in  wind  velocity  can 
be  examined  with  the  aid  of  Equation  (17).  The  use  of  this  equation  is  quite  relevant  since 
the  shell  spends  a  large  amount  of  time  near  apogee,  and  thus,  the  trajectory  should  be 
significantly  afFected  by  changes  in  atmospheric  values  near  apogee.  The  variation  of  its 
natural  log,  when  set  to  zero,  is  equivalent  to,  for  instance,  changing  the  density  and  wind 
velocities  so  that  the  acceleration  of  the  projectile  remains  unchanged.  When  the  down-range 
wind  velocity  and  a  are  free  to  vary,  the  resultant  trajectory  will  not  vary  from  its  original 
path  if 


6a  26Wi 
~a  ~  Ui  -  Wi  ‘ 


(22) 
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6.  TEST  RESULTS  AND  DISCUSSION 


The  trajectory  was  fitted  along  its  intervals  in  a  sequential  fashion,  starting  at  ground 
level  and  using  the  ending  values  of  the  intervals  as  an  initial  estimate  for  the  next  interval. 
Although  the  initial  estimates  also  serve  as  the  unchanging  initial  values  for  the  density  and 
temperature,  for  the  wind  velocity  values,  initial  estimates  are  adjusted  to  obtain  the  best 
fit  along  with  the  fitted  lapse  rate.  The  wind  velocity  data,  together  with  the  curves  found 
by  the  least  squares  fitting  method,  are  shown  in  Figure  4. 


Figure  4.  Comparison  of  Range  Wind  Velocity  (Wj)  and  Cross  Wind  Velocity  (W3)  Data 
with  Respective  Fitted  Values 


The  profiles  are  roughly  reminiscent  of  boundary  layer  behavior,  as  one  would  expect.  The 
wind  velocities  are  low  near  the  ground  and,  at  first,  generally  increase  with  altitude.  The 
agreement  of  the  least  squares  fitting  results  with  the  data  is  good,  especially  considering 
that  a  constant  value  of  wind  velocity  is  assumed  for  each  segment.  These  fitted  values, 
although  not  fitted  on  the  same  intervals  as  are  used  for  firing  tables,  are  in  a  format  similar 
to  the  standard  MET  message. 
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The  density  and  temperature  data  are  compared  with  the  fitted  data  in  Figure  5.  The 
fitted  density  agrees  best  with  the  data  near  the  ground  and  agrees  less  well  near  apogee 
since  the  calculation  is  started  at  the  lowest  lying  segment  and  proceeds  stepwise  upward 
with  the  ending  density  value  for  the  lower  segment  becoming  the  initial  value  for  the  higher 
segment.  Errors  made  near  the  beginning  can  propagate  and  amplify  over  the  segments  at 
higher  altitudes.  In  fact,  the  fitted  density  values  appear  to  diverge  from  the  data  at  the 
higher  altitudes. 


ALTITUDE  (M) 

Figure  5.  Comparison  of  Density  and  Temperature  Data  with  Fitted  Values 

The  temperature  values  obtained  throughout  the  fitting  calculations  are  erroneous  by 
a  few  degrees  over  most  of  the  trajectory.  The  initial  value  for  the  temperature  is  taken 
on  the  ground.  A  temperature  inversion,  which  often  occurs  in  the  winter,  results  in  a 
rapid  initial  increase  in  error  with  height.  An  error  in  the  temperature  chiefly  affects  the 
value  of  the  zero-yaw  drag  coefficient,  Cdq,  which  depends  weakly  on  the  temperature  for 
supersonic  velocities,  more  strongly  for  transonic  velocities,  and  not  at  all  for  the  low  to 
moderate  subsonic  velocities.  A  major  part  of  the  trajectory  occurs  at  subsonic  velocities. 
Nevertheless,  the  first  part  of  the  trajectory  occurs  at  supersonic  and  transonic  velocities. 
The  error  in  the  temperature  values  also  affects  the  density  values  since  the  change  in  the 
density  depends  upon  the  change  in  the  temperature. 
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ERRORS  INDUCED  BY  MEASUREMENT  ERRORS  AND  MODEL 
ASSUMPTIONS 


With  the  assumptions  of  a  constant  wind  velocity  and  a  linear  temperature  variation 
within  a  segment,  it  is  obvious  that  the  fitted  trajectory  could  only  approximate  the  trajec¬ 
tory  that  was  generated  by  using  real  atmospheric  data.  Figure  6  shows  the  positional  errors 
in  meters  as  a  function  of  altitude.  The  positional  error  is  simply  the  difference  between 
the  fitted  value  of  the  position  and  the  generated  value  of  the  position.  Here,  a  large  error 
occurs  at  a  lower  altitude  with  no  spiking  thereafter  except  for  the  maximum  errors  at  the 
segments  that  occur  at  their  end  points.  This  large  enor  at  the  lower  altitudes  is  probably 
associated  with  the  projectile  being  in  the  supersonic  and  transonic  regions.  If  an  error  is 
made  in  obtaining  the  temperature,  the  drag  coefficient  can  be  significantly  affected  while  in 
the  subsonic  region  where  the  drag  coefficient  is  only  weakly  dependent,  on  the  Mach  number. 
Errors  near  the  beginning  of  the  segment  must  be  small;  otherwise,  errors  made  near  the 
beginning  of  the  trajectory  over  the  segment  would  tend  to  propagate  and  increase  over  the 
sc^rncnt . 

°  XI  DIFFERENCE 

X2  DIFFERENCE 


As  discussed  earlier,  there  will  be  an  uncertainty  in  the  measurements  of  the  end  posi 
tions  of  a  trajectory  segment  which  will  induce  errors  in  obtaining  the  fitted  values  of  the 
meteorological  quantities.  Uncertainty  in  measured  positions  along  the  trajectory  segment 
should  also  induce  errors  in  fitting  the  atmospheric  values.  These  uncertainties  are  explored 
here  by  superimposing  Gaussian  noise  with  a  standard  deviation  (SD)  of  0.2  meter  distance 
upon  the  trajectory  position  points.  To  obtain  a  more  detailed  view  of  the  error,  we  turn 
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to  numerical  analysis.  It  is  .shown  that  the  number  ot  significant  digits,  J,  in  the  density  is 
given  as 

j  =  log.o  ,  (23) 

in  which  Ep  is  the  relative  error  of  the  density.  For  instance,  Equation  (23)  shows  that  if 
the  relative  error,  Ep,  has  the  value  0.05,  the  nrunber  of  significant  digits,  J,  is  equal  to  one. 
As  the  relative  error  becomes  smaller,  the  number  of  significant  digits  becomes  larger.  The 
number  of  significant  digits  of  resultant  relative  error  in  density  for  the  fitted  values,  for 
both  noise  and  no  noise,  are  shown  in  Figure  7.  Convergence  to  a  solution  took  much  longer 
with  the  noise  present.  Not  surprisingly,  the  fitted  solution  with  no  nobe  results  in  a  better 
fit.  It  also  appears  that  the  errors  obtained  at  lower  altitudes  may  propagate  to  degrade 
the  solution  at  higher  altitudes.  As  discussed  earlier,  however,  longer  segments  are  needed 
)  obtain  accurate  fits  as  the  altitude  increases.  Fitting  the  trajectory  with  Gaussian  noise 
with  a  SD  distance  of  1  meter  resulted  in  poor  agreement  and  is  not  shown  here. 
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Figure  7.  Significant  Digits  of  Relative  Error  in  Density  vs.  Altitude. 

Another  source  of  uncertainty  is  the  actual  value  of  the  drag  coefficient,  C p.  The  value 
of  the  drag  coefficient  is  commonly  known  to  within  only  a  half  percent,  as  discussed  earlier. 
The  drag  coefficient  value  may  also  vary  even  more  from  lot  to  lot.  To  test  this  source 
of  uncertainty,  the  fitted  value  of  the  drag  coefficient  was  assumed  to  be  a  half  percent 
lower  than  the  drag  coefficient  values  used  in  generating  the  trajectory.  The  resulting  ntted 
density  was  approximately  a  half  percent  higher  than  had  been  obtained  previously.  This 
result  would  be  expected. 

The  approach  of  using  normally  stale  density  and  temperature  data  to  obtain  equivalent 
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winds  can  also  be  briefly  scrutinized  and  perhaps  be  put  in  perspective.  Table  2  shows  the 
SD  of  atmospheric  quantities  measured  at  selected  times  after  an  initial  measurement.  This 


Table  2.  Errors  Incurred  in  Stale  Met  Messages  -  M483A1  Shell 


Staleness 

&wind 

G den 

&  lemp 

^ron 

°dejl 

Hours 

m/s 

% 

% 

m 

m 

0 

r'0:4~ 

0.15 

0.25 

14 

8 

1 

2.1 

0.40 

0.30 

66 

40 

2 

2.5 

0.69 

0.57 

84 

49 

4 

0.97 

0.79 

122 

71 

oo 

m 

6.60 

3.0 

341 

109 

table  was  generated  from  a  report  by  Matts  and  F  is  (1990).  Shown  also  are  the  associated 
error  magnitudes  for  range  and  deflection  for  the  M483A1  shell  fired  from  the  M109A2 
howitzer.  The  target  range  is  15,000  meters  with  the  zone  8  charge. 

This  table  shows  nonzero  values  with  fresh  atmospheric  data.  These  values  are  equal  to 
the  instrument-measuring  errors.  The  errors  incurred  by  not  taking  into  account  atmospheric 
conditions  are  large,  as  shown  in  the  last  row  of  Table  2.  If  we  assume  that  the  density  and 
ternpe-ature  could  be  approximated  by  a  standard  atmospheric  model,  then  using  a  shell’s 
measured  trajectory  to  find  the  wind  velocity  by  a  fitting  process  would  yield  huge  fictitious 
wind  velocities.  If  we  use  density  and  temperature  data  that  are  a  few  hours  old,  however, 
we  can  use  Equation  (22)  to  find  the  projectile  velocity  necessary  to  leave  the  net  force  of 
the  projectile  unchanged.  For  a  staleness  time  of  4  hours,  the  value  of  the  shell  velocity  to 
satisfy  Equation  (22)  is  approximately  740  m/s.  This  velocity  will  be  higher  than  the  the 
actual  velocity  of  the  shell  at  apogee  and  shows  that  the  error  in  the  density  contributes 
more  to  the  range  error  than  is  contributed  by  the  error  in  wind  velocity.  Pedersen  (1994) 
is  exploring  the  possibility  of  using  more  recent  ground  observations  to  make  an  assertion 
about  how  atmospheric  conditions  have  changed  at  higher  altitudes.  Such  knowledge  would 
avoid  the  large  errors  in  range  and  deflection  that  might  occur  if  1 1  >  artillery  piece  were  fired 
on  a  different  azimuth. 

Measurements  over  the  entire  trajectory  have  not  been  available,  and  particularly,  mea¬ 
surements  have  not  been  available  over  the  trajectory  portion  near  the  ground.  Although  the 
current  measurements  appear  as  a  smooth  curve,  the  smoothness  is  an  artifice  of  the  data 
reduction  technique.  For  both  CPS  and  radar  systems,  the  actual  uncertainty  can  be  meters. 
According  to  the  analysis  made  earlier,  the  uncertainty  in  the  atmospheric  quantities  would 
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be  unacceptably  large  for  the  segment  lengths  used  here. 

Other  approaches  are  obvious  to  obtain  sufficiently  accurate  atmospheric  values  but  will 
only  be  discussed  here.  Accuracy  should  be  improved  by  using  the  whole  trajectory  length 
although  the  calculational  effort  would  be  much  increased.  The  atmospheric  values  would 
have  to  match  at  any  given  altitude  while  both  ascending  and  descending,  thus  increasing 
the  number  of  constraints  on  the  solution.  The  density  should  be  fitted  as  a  continuous 
function  while  the  fitted  wind  velocities  would  be  assumed  as  constant  values  over  segments 
of  the  trajectory.  This  method  would  be  slow  to  converge  to  a  solution,  and  the  accuracy 
of  the  fitted  values  might  still  be  affected  by  the  drag  coefficient  and  the  velocity  appearing 
in  the  flight  dynamic  equations  as  complementary  factors.  As  for  the  current  approach,  a 
bad  guess  for  the  initial  estimates  may  result  in  a  converged  solution  that  is  far  from  the 
absolute  minimum. 

8.  SUMMARY,  CONCLUSIONS,  AND  RECOMMENDATIONS 

By  using  a  least  squares  fitting  procedure,  we  obtained  atmospheric  conditions  knowing 
only  the  projectile’s  trajectory  and  its  aerodynamic  coefficients.  The  projectile  trajectory  was 
generated  using  the  MPM  model.  The  atmospheric  conditions  were  supplied  by  a  weather 
balloon.  A  least  squares  fit  of  several  trajectory  segments  that  comprised  altitude  zones  was 
made.  The  fitting  procedure  started  from  the  ground  and  used  the  fitted  density  and  temper¬ 
ature  conditions  at  the  end  point  of  the  trajectory  segment  as  the  initial  fitted  conditions  for 
the  next  segment.  The  agreement  of  the  fitted  meteorological  values  with  the  balloon  data 
was  good,  but  the  agreement  decreased  with  the  increasing  amounts  of  noise  added  to  the 
trajectory  positions.  A  superimposed  noise  with  a  SD  of  1  meter  resulted  in  an  unsatisfac¬ 
tory  fit.  These  results  agree  with  an  error  analysis.  To  obtain  atmospheric  conditions  with 
the  needed  precision,  the  accuracy  study  also  shows  that  the  starting  and  ending  positions 
of  segment  lengths  must  be  known  with  more  precision  than  current  measuring  techniques 
are  capable  of  providing.  A  different  approach  is  needed,  such  as  using  the  data  over  the 
whole  trajectory  instead  of  fitting  to  the  data  from  a  relatively  small  trajectory  segment. 

The  calculaiional  effort  would  be  dramatically  reduced  if  the  temperature  were  known, 
and  hence  the  temperature  dependent  flight  coefficients,  such  as  the  drag  coefficient,  would 
be  immediately  known.  A  thermometric  device  on  board  a  flight  projectile  could  radio  the 
temperature  data  to  a  ground  station.  Although  the  hardware  complexity  would  be  some¬ 
what  increased,  this  approach  might  yield  a  practical  robust  solution  with  rapid  convergence. 
Calculations  could  be  accelerated  even  more  if  the  density  could  be  determined  from  the  as¬ 
sumption  of  a  hydrostatic  atmosphere,  knowledge  of  the  temperature,  and  the  equation  of 
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state.  To  obtain  the  equation  of  state,  however,  the  moisture  content  of  the  air  must  be 
known.  Generally,  the  moisture  content  of  the  air  is  appreciable  near  the  ground.  Perhaps, 
one  could  use  the  current  moisture  content  on  the  ground  together  with  older  detailed  data 
to  supply  data  with  sufficient  accuracy. 
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LIST  OF  SYMBOLS 


a  AKV/t) 

A  (£fc»)  (t)  •  M 

AZ  azimuth  of  the  1-axis,  measured  clockwise  from  North 


(£)©’• 


c 

Coriolis  acceleration,  — 2u;  x  U 

CD 

drag  coefficient:  | drag  force |  =  (pV2S/2)Co 

Cooi  C-D2 

zero-yaw  and  yaw-drag  coefficients: 

Cd  —  Coo  +  Cd2\&A2 

Ct, 

roll  damping  moment  coefficient: 

| roll  damping  moment\  =  ±(pV2S(/2)(<j>£/V)Ce] 

cLa 

lift  force  coefficient: 

1  Hfi  force \  =  ±{pV2S/2)\ae\CLa 

Cm „ 

static  moment  coefficient: 

\static  moment\  =  ±(pV2S£/2)\ae\CMa 

Cnp„ 

Magnus  force  coefficient: 

\Magnus  force]  =  ±(pV2S/2){<j>e/V)\ae\CNpa 

CU...C3  fitting  parameters  defining  the  wind,  Eqs.(9)  and  (11) 

D  [l/»] 

E  gun  elevation 
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relative  error  in  the  fitted  air  density 


gravity  acceleration 


at  sea  level;  the  STAT  value  at  45  deg.  N.  latitude  is  9.80665  m/s 


g  -f  C,  gravity  plus  Coriolis  acceleration 


_JL_  \q  + 

1  Wa  ^  » 


[m/s2] 


G i,  G2,  G3  flat-earth  system  components  of  G 


(P) 

v  (§&)  w 


geopotential  altitude,  Eq.(12) 


value  of  H  at  the  starting  time  of  a  fitting  interval 


value  of  H  at  the  bottom  of  a  STAT  altitude  zone  (Table  1) 


axial  moment  of  inertia 


number  of  significant  digits  in  the  fitted  air  density 


Ix/rn(2 


reference  length 


latitude  at  the  launch  point  (for  Southern  Hemisphere 
firings,  replace  L  by  -L) 
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m 

M 

Q 

R 

n 

s 

STAT 

t 

i 


r 


v 

TA 

Tb 

u 


Uu  U2,  lh 


V 

V 

V. . 


projectile  mass 
molecular  weight 

constant  in  the  STAT  air  density  formula,  (13) 
effective  radius  of  the  earth  (6,356,766  m) 
universal  gas  constant 
reference  area,  iti2 /4 

1966  U.S.  Standard  Atmosphere,  referenced 
time 

nondimensionalized  time,  U\b{t  —  h)/t 
temperature,  degrees  Kelvin 

temperature  gradient  dT/dH;  fitting  parameter  C3 

temperature  at  the  starting  time  of  a  fitting  interval 

STAT  temperature  at  altitude  Hb 

projectile  velocity  with  respect  to  the  earth 

**♦ 

flat-earth  system  components  of  U 

\v\ 

U  —  W ,  projectile  velocity  with  respect  to  the  air 
speed  of  sound 
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Vi,  v2,  v3 


fiat-earth  system  components  of  K 


W  wind  velocity  with  respect  to  the  earth 

W\,  W2,  W3  flat-earth  system  components  of  W 

X  projectile  position  with  respect  to  the  earth 

x  nondimensionalized  position  variable,  (A”i  —  Xn)/i 

X\,  X2,  X3  flat-earth  system  components  of  X 

X\  :  down  range 

X2  '■  the  height  above  sea  level 

X3  :  lateral,  positive  to  the  right  looking  down-range 


0 

P 

Pa 

Pb 

4> 

u 

U>e 


,  W'i,  U>3 


(') 


the  yaw  of  repose 

the  mass  fraction  of  water  vapor  in  the  air 
air  density 

air  density  at  the  starting  time  of  a  fitting  interval 

STAT  air  density  at  Hb ,  see  Table  1 

axial  spin  rate 

angular  velocity  of  the  earth 

H 

flat-earth  components  of  u5: 

ue{cos  L  cos  AZ,  sin  L,  —  cos  L  sin  AZ) 


<n  )/* 
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